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Particles added to a fluid interface can be used as a surface stabilizer in the food, oil and cosmetic 
industries. As an alternative to rigid particles, it is promising to consider highly deformable particles 
that can adapt their conformation at the interface. In this study we compute the shapes of soft 
elastic particles using molecular dynamics simulations of a cross-linked polymer gel, complemented 
by continuum calculations based on linear elasticity. It is shown that the particle shape is not only 
affected by the Young’s modulus of the particle, but also strongly depends on whether the gel is 
partially or completely wetting the fluid interface. We find that the molecular simulations for the 
partially wetting case are very accurately described by the continuum theory. By contrast, when 
the gel is completely wetting the fluid interface the linear theory breaks down and we reveal that 
molecular details have a strong influence on the equilibrium shape. 


I. INTRODUCTION 


An important application of particle stabilized fluid interfaces goes back to the importance of retaining 
the dispersivity of emulsions 1 ®. The dispersion process results in a large interfacial area and hence high 
interfacial energy. Being energetically unfavourable, dispersed phases will eventually coarsen and form phase 
separated volumes of fluidffl A traditional way to stabilise the dispersions with respect to phase separation 
is to use surfactants. As an alternative for surfactants, an interface can also be made kinetically stable 
by adding solid particle^ 1 ®. Particle-stabilised emulsions, also called Pickering emulsions® are metastable 
since the particles anchor to the interfaces much more strongly than the surfactant molecules. For example, 
particles as small as a few tens of nanometres have a desorption energy as high as 10 3 to 10 4 fcsT, while it 
is around 10 ksT for surfactant molecule^®. 

As an alternative to solid particles, soft particles have recently attained attention as stabilizers for 
emulsioniP. The shape of the particles can adapt to the interface and depends on the interplay of the 
molecular interactions between particle and fluids as well as the elastic properties of the particle itself. The 
most important types of releva nt par ticles are cross-linked polymer networks (microgels(pH 12 * and l igand 
or polymer grafted nanoparticles 13 15; . Combinations of both types are so-called core-shell microgels 16 17 . 
The desorption energy of micron-sized microgel particles at an interface was found to increase even up to 
10 6 A;.bT®, making them very promising as stabilizer 

The presence of a network of interconnected polymer chains makes these particles deformable under an 
external forcing. They can show the properties of both the individual chains as well as the particles with well- 
defined boundaries®. This highlights the role of the interfacial tension of soft particles in their deformation 
at a fluid interface. Understanding the details of the adsorption of soft particles to fluid interfaces is highly 
complex, as it involves e.g. the cross-linking of the polymer network, temperature and pH value 22 . However, 
obtaining a fundamental understanding of the conformation of these objects at a fluid interface in terms of 
the macroscopic and microscopic parameters is crucial for exploiting them in many practical applications. 

Most experimental studies on the interaction of so ft particles and fluid interfaces focus on the behaviour 
of a single microgel particle at a fluid interface®®*! or bulk properties of microgel covered interface^® 31 
whose key observation is that the microgel particles take a core-corona (also called fried-egg) shape at a fluid 
interface. There are some challenges when studying the deformation of a particle at a fluid interface: due 
to the small size of microgel particles and their reflective index being very close to that of the surrounding 
solvent, imaging is difficult®. In addition, soft particles are very sensitive to external stimuli such as 
pH or temperature change. There are few numerical simulations of polymer grafted nanoparticles using 
molecular dynamics. For example, Schwenke et al. 32 have studied the particle conformation and particle- 
particle interaction of polymer coated nanoparticle in different solvents. Udayana et al!® have studied 
the adsorption/desorption energy of polymer coated nanoparticles to/from an interface. Lane et al!® have 
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liquid droplet soft particle 


FIG. 1 Snapshots from molecular dynamics simulations of soft particles at a fluid-fluid interface. The key 
parameters of this study are the elastic modulus and the wettability of the particles. The left panels show the 
outcome for a liquid drop, consisting of a non-cross-linked polymer, respectively for partially wetting (top) and 
complete wetting (bottom). The right panels show the corresponding particle shapes when the polymers are 
cross-linked and have a finite elastic modulus. The reason for small asymmetry between left and right hand side of 


the soft particle in the complete wetting regime is explained in section IV B 


investigated the distribution of polymer chains. None of these works has discussed the shape of the soft 
particles in terms of their macroscopic properties. 

From a theoretical point of view, it is appealing to investigate the shape of soft particles using a continuum 
framework based on linear elasticity. As the particle becomes increasingly soft, however, such a theory must 
also account for interfacial forces. The interplay of capillarity and elasticity was recently investigated in 
great detail in the context of adhesion of particles®^® or the wetting of liquid droplets on highly deformable 
substrateIn these studies the governing parameter was found to be the elasto-capillary length 7 a /E, 
comparing the surface tension of the solid % to its Young’s modulus E. If the particle (or droplet) size 
R is large with respect to 7 a /E, it can be considered effectively rigid - except in a small region near the 
contact line. In the opposite limit where the particle is small with respect to the elasto-capillary length, the 
elasticity is so weak that the elastic medium can effectively be considered as a liquid with a surface tension. 
Extending this point of view to the adsorption of soft particles at interfaces, one thus expects two limiting 
cased®® one is the behaviour of a perfectly rigid particle at an interface (governed by the Young contact 
angle), while the soft extreme of vanishing elasticity corresponds to a liquid droplet at a fluid-fluid interface 
(governed by the Neumann contact angles 46 ). It has remained unclear to what extent linear elasticity can 
describe the shape of soft particles at fluid-fluid interfaces, in particular when the surface tension of the 
particle is relatively low®. 

In this paper we quantify the equilibrium shapes of soft particles at fluid-fluid interfaces, by combining 
molecular dynamics simulations and exact solutions derived from linear elasticity. Previously, the molecular 
dynamics method has been used to study the interaction between a liquid droplet or a solid nanoparticle 
with a fluid interface®' , a nd the inte raction between soft materials made of inter-connected polymer chains, 
with liquid or solid surfaced®®®® Similarly, our particles consist of a cross-linked polymeric liquid and 
are adsorbed at a fluid-fluid interface. By varying the molecular interactions and cross-linking density, we 
can explore a broad range of Young’s moduli and interfacial energies. Snapshots of typical simulations are 
shown in Fig. [TJ where the individual polymer chains can be identified by their color. Apart from an expected 
dependence on the Young’s modulus, governed by the dimensionless parameter 7 a /(ER), our key finding 
is that one needs to distinguish the cases where the polymer is partially or completely wetting. The left 
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column of Fig. [ljshows simulations of the polymer liquid (i.e. without cross-linking), respectively for partial 
wetting (upper panel) and complete wetting (lower panel). The right column shows equivalent systems, but 
now for soft solid particles that are cross-linked and thus exhibit a finite elastic modulus. We find that 
the partially wetting systems are very accurately described by linear elasticity, as long as the solid surface 
tension is comparable or larger than that of the fluid-fluid interface. In the case of complete wetting, by 
contrast, the interfacial forces favor a state where the gel covers the entire interface but this is prohibited 
by the network elasticity. This leads to an intricate elasto-capillary balance, where the elastic response is 
highly nonlinear. In addition, we observe an important influence of molecular details of the cross-linking, 
signalling a breakdown of continuum theory. 

The paper is organised as follows: First, we discuss the analytical continuum theory, which is limited to 
the regime of small deformations, in section [TTJ Details of the molecular dynamics simulations are presented 
in section [Til] and our simulation results are presented in section [W] Here we provide a detailed comparison 
to the continuum predictions and explore the regime of complete wetting. The paper closes with a discussion 
where we also comment on the implications of our findings for microgel particles (section |V| . 



FIG. 2 Schematic shape of the undeformed (grey circle) and deformed particle (solid line), (r, <f) denotes the 
components of the polar coordinates. The fluid-fluid interface is located along the x-axis, at (f> = 0 and n. R, A R 
and 6 S denote the particle’s radius, its maximum deformation and contact angle. 


II. MACROSCOPIC THERMODYNAMIC FORMULATION 


In this section we pose the problem of an elastic particle at an interface from a macroscopic thermodynamic 
perspective. The relevant dimensionless parameters are identified and the shape of the deformed particle is 
computed in the framework of linear elasticity. This provides a benchmark for the case of weakly deformed 
particles, and outlines the key issues that will later be explored by molecular dynamics for large particle 
deformations. This study is done in a two-dimensional setup. It is expected that the two-dimensional 
simulations capture the salient features of the three-dimensional configuration because the main physics 
of the problem happens at the vicinity of the contact lines except that the global deformation in three- 
dimensions will be slightly less than in the two-dimensional case due stress in the azimuthal direction. 


A. Dimensionless parameters 

We consider a two-dimensional elastic particle (Young’s modulus E) that is placed at a liquid-liquid 
interface (surface tension 7^), as sketched in FigJT] For simplicity, the study is limited to the cases where the 
solid particle has an equal surface tension with the two liquid phases, simply denoted 7 S , so that the particle 
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will develop a symmetric shape 


Notel 


The reference state of the particle is a circle of radius R, which is for 


example achieved by cross linking a liquid drop that is fully immersed in one of the two liquids. We assume 
the particle to be incompressible (Poisson ratio v = 1/2), although the results are easily general ised to an 
arbitrary Poisson ratio. In analogy to the case of liquid drops on deformable substrate d 40 * 41 ! 52 !, there are 
two dimensionless parameters that characterise the problem, namely 


ci 7s , 7s 

S'=—, and 7 =—. 
ER 7f, 


(1) 


The first term is called the “softness” of the particle and compares the elasto-capillary length ^f s /E (based on 
the solid surface tension), to the particle radius. The second term is the ratio of surface tensions, controlling 
the wetting conditions. Before turning to a detailed analysis, it is instructive to discuss the extreme limits 
of the softness S. For S = 0, the particle can be considered as undeformable (large E or R) and one recovers 
the usual case of a rigid spherical colloid at an interface. The limit S = oo corresponds to a particle without 
any rigidity, as is achieved for a polymeric drop without any cross-linking. Hence, in this limit one expects 
to recover the “liquid lens” floating at the interfac^^. The shape of such a lens is governed by the ratio of 
surface tensions 7 , and its contact angle 6 S follows from the Neumann triangle (Figj2]). Given the symmetry 
of our problem, where the gel has an equal surface tension with the two liquids, the Neumann balance of 
surface tensions reads 


7 b = 27 s cos 9 S 


9 S = arccos 



Quantifying the maximum drop deformation A R as defined in Fig(2j we find that 


( 2 ) 


A R, 


•drop 


R 


29 s — sin 29. 


■ sin 9 s — 1 for S 


00 . 


( 3 ) 


Note that such a drop exhibits a transition from partial wetting to complete wetting when the surface 
tension ratio 7 —> 1/2: at smaller values of 7 the Neumann balance ([ 2 ]) cannot be satisfied. The transition 
corresponds to 9 S —> 0 for which Eq. ([3]) gives a divergence of A R/R, signalling the onset of a wetting layer. 


B. Small deformations: linear elasticity 

In the limit of small deformations, one can resolve the particle shape using linear elasticity theory. The 
two-dimensional problem consists of a disk of radius R that is deformed under the influence of an external 
stress at the free surface. This stress, or traction, has two contributions. First, the liquid-liquid surface 
tension pulls on the particle at <f> = 0 and n, where we use polar coordinates as defined in Fig. [2] This 
traction is described by a perfectly localised force per unit length, using a Dirac <5 distribution. Second, the 
deformation changes the curvature of the particle-liquid interface, and as such induces a capillary pressure 
that acts as an additional traction r y s n. In this expression is the surface stress (here simply denoted as 
surface tension), while n is the extra curvature of the solid-liquid interface due to the deformation, which to 
linear order in the radial displacements u r reads 

1 ( d 2 u r \ 

K = w{ Ur + w)' (4) 

to be evaluated at the disk boundary r = R. The importance of this solid Laplace pressure has been 
highlighted in several recent papers involving very soft geld^^, and will also be apparent in the present 
work. The total traction boundary condition thus becomes 


Notel j n case the surface energy varies with the amount of stretching, one needs to distinguish between surface free energy and surface 
stress for the liquid-gel interface due to the Shuttleworth effeclP^l Since here we consider a system with two identical liquids, 
this effect does not induce a forcing tangential to the particle surfacd^ and our calculation is valid - in the remainder we 
simply refer to surface tension, which for the liquid-gel interface must be seen as the surface stress. 
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<r rr (r = R,(/>) = ^ [£((/>) + 5{(j) ~ ?r)] +7 s«, (5) 

cr r <t>{r = R, $) = 0 , ( 6 ) 

where the latter expresses the no-shear stress boundary condition. Interestingly, the traction depends on 
the displacement. As a consequence, the resulting tractions and displacements have to be determined self- 
consistently. 
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FIG. 3 Deformation of a 2 D particle at an interface, predicted by equation (|tJ) for 7 = 2.0 and S = 200000, 0.2 and 
0.02. The contact angle 8 S = 75.5° for all cases, dictated by the Neumann balance. 


This two-dimensional problem can be solved using the Airy stress function forrnalisnP23, and the complete 
derivation is given in the Appendix. The result is obtained using a series expansion in polar coordinates, 
based on the classical solutions by Michell. The resulting shape of the deformed particle is given by (cf. 
Appendix) 


R 

which indeed depends on both the dimensionless softness S and the ratio of surface tensions 7 . A similar 
expression was derived for axisymmetric particles in 45 . Some examples of the deformed particles are shown 
in Fig. [3] Importantly, the particle exhibits a finite contact angle d s , for all values of the softness S, as 
can be seen from the kinks at the position where the fluid-fluid interface applies its force. In the theory, 
the angle 9 is even independent of the elastic properties of the particle. This can be inferred from the large 
k -asymptotics of ([T]), which gives a slope discontinuity in dr p /dcj) of magnitude I/ 7 . In terms of the contact 
angle this gives 


i + — V 


4 k cos 2 k<f> 


2'yn jy-J (4 k 2 — 1)(1 + 3 Sk) ’ 


( 7 ) 


9s 


7 r 1 

2 _ 27 ’ 


( 8 ) 


which is the small slope, or large 7 , equivalent of (J2|. Hence, the continuum theory predicts that the contact 
angle is determined solely by the surface tensions, satisfying Neumann’s law, regardless of the stiffness. This 
is in perfect analogy to results previously obtained for liquid drops on soft solids. The Neumanns angle also 
appear for axisymmetric soft particles 45 . 

The maximum radius of the particle is attained at <f> = 0 and 7 r, for which the series |7| can be summed 
to the explicit form 
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A R _ 35 [2 - 2'k(l/35) - 41n(2) - 2 lE - 35] 

~ 7(95 2 - 4)tt ' 

Here 7 ^ is the Euler-Mascheroni constant and 4' is the digamma function. As expected, the deformation 
A R/R vanishes for a rigid particle (5 = 0) and saturates to a finite value in the very soft limit (5 —► 00 ). 
In the soft limit, we find A R/R —► 1 /(^ 7 ), which agrees with 0 for large 7 . It the rigid limit, we obtain 
[using ^(a;) ~ In a;]: 


A R 3 5, /a\ 

^~2iy ,n \s) 


4 7 b , 

Her'" 


( 7 ) 


for S<1. 


( 10 ) 


where the numerical constant a = |e 7js_1 . The deflection thus increases linearly with softness, involving 
the fluid-fluid surface tension 75 , with logarithmic corrections that involve the the su rface ten sion of the 
solid 7 s. Similar scaling laws for contact line deflections were obtained for drops on gel ci 40 * 44 ^ . Hence, the 
surface tension of the solid 7 S is critical to achieve a finite elastic deformation in the macroscopic theory: 
otherwise, the deformation would be (logarithmically) divergent. In addition, 7 S determines the boundary 
condition for the contact angle, following Neumann’s law. We will comment later how this picture must 
breakdown in the regime of complete wetting, for which the surface tensions cannot achieve the Neumann 
balance. 

The results from the molecular dynamics simulations will be directly compared to the shape 0 , the 
Neumann contact angle prediction 0, as well as the maximum extension A R. For the latter, it is convenient 
to normalise the extension § by its limiting value achieved for a liquid droplet at large 7 i.e. 1 /(^ 7 ), so 
the result only depends on the softness 5 (not on 7 ): 


A R _ 35 [2 - 2T(l/35) - 41n(2) - 2j e - 35] 
A Rdrop ~ (95 2 - 4) 


( 11 ) 


C. Observations and questions 

The macroscopic theory provides quantitative predictions, for the shape and the contact angle, that will be 
tested using molecular dynamics simulations. However, the simulations are even of more interest when the 
linear theory breaks down, and we therefore point out some important limitations of linear elasticity. First, 
the analysis only applies for large 7 , for which according to 0 the Neumann angle 9 S is always close to 7 t/ 2 . 
The reason for this is that in the linear elasticity calculation we have evaluated the curvature k to lowest 
order in du r /d(f> - keeping nonlinear terms prevents the solution by Fourier expansion. As a consequence, 
the theory cannot capture the wetting transition at 7 = 1/2, for which the angles are no longer small. Of 
particular interest is the force balance in the vicinity of the contact line. Contrarily to the linear theory, 
for the situation of complete wetting the Neumann balance cannot be satisfied and an elastic contribution 
must emerge to balance the localised contact line force. How does the wetting transition appear in the case 
of finite elasticity? 


III. MOLECULAR DYNAMICS SIMULATION 

We study the particles in a quasi two-dimensional setup, allowing to simulate a large particle at a less 
computational cost as well as a direct comparison to the elasticity theory described above. The setup is 
prepared in three main steps. First a fluid-fluid interface is made. Then a cylindrical particle is built from 
polymeric chains and finally the particle is inserted at the interface and the system is equilibrated. In 
addition, two auxiliary setups are needed to measure the Young’s modulus of the particle and the liquid- 
particle surface tensions. Below, these steps will be outlined in detail. We use the molecular dynamics 
method implemented in the GROMACS 5.0. software to do the numerical simulations. Integration of 

Newton’s equations is done using the leap-frog algorithm. Neighbour searching is done using verlet dynamic 
lists. Visualization of the results is done using the VMD packag^l 
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A. Interface modelling 


To keep the model as simple as possible, we use the modified shifted and truncated Lennard-Jones 
potential to model the interaction between each pair of particles, 


U = 



if t < r c 
if r > r c . 


( 12 ) 


Here r is the distance between two particles, d is the repulsive core diameter fixed to 0.34 nm, e is the depth 
of the potential minimum set to 3 kJ/mol, and r c is the the cut-off radius equal to 5 d. The characteristic 

time of the system is r = yj' r ^L = o. 6 ps and m = 10 amu is the particle mass. Using this set of parameters, 
the surface tension of the liquid-liquid interface would fall in the range between 0 and 73 mN/m, which are 
typical experimental values. In the remainder of the paper, we will not use SI units, but reside to either plots 
nondimensionalized by the particle radius, or use r, d, and e as the reference values for the time, length, 
and energy, respectively. The coefficient a is set to 1 for two identical Lennard-Jones particles and it is 
smaller than 1 when two particles are of different species. This naturally leads to phase separation, where 
a controls the degree of miscibility and hence the surface tension. In this study there are three phases and 
hence three values for a are required. Since only symmetric liquid-particle interactions are considered, both 
fluids have the same interaction parameter with the gel, which we call a p . Therefore, only two values of a 
and a p will be varied to produce different values for the ratio of surface tensions 7 . 

The liquid-liquid interface is produced under the NP n T ensemble where P n denotes the pressure normal 
to the interface^. For keeping the temperature and normal pressure constant, Nose-Hoover temperature 
coupling^ and Parrinello-Rahman pressure coupling^ methods are used, respectively. Periodic boundary 
conditions are applied in all three directions. In our simulations, the produced liquid-liquid interface has 
approximately two million Lennard-Jones particles and its dimensions in x , y, and z directions, as defined 
in Fig. [TJ are 600d, 300d, and 14.7 d, respectively. 


B. Particle modelling 

The particle is made from interconnected polymeric chains. Each chain is constructed using a coarse¬ 
grained model of beads and springs and consists of 32 monomers where two neighbouring beads interact via 
the following potential: 

U{r) =-\k s R 2 max —) + 4 e (r) . (13) 

The first term on the right hand side is the so-called FENE (finite extensible nonlinear elastic) potential 
which controls the amount of the extension. The second term is the Lennard-Jones repulsion term that 
accounts for the reduced volume effect and prevents the collapse of the beads onto each other. We fix 
k s = 25 (^ 2 ) and R m a X = 1-5<L These values are chosen such that the simulation timesteps can be taken 
on the order of the pure Lennard-Jones fluid while making the link breakage energetically impossible 67 . 
By randomly cross-linking the polymeric chain, one creates a network of entangled polymers that exhibits 
elasticity in the long time limit. By increasing the density of the cross-linking, the rigidity of the network 
increases leading to a larger Young’s modulus. All non-neighbouring beads of the chains interact with each 
other and also with the beads of the other chains through the Lennard-Jones potentials as described in the 
previous section. It is assumed that the particle is an isotropic material and the density of the cross-links 
is distributed uniformly over the particle volume. In addition, the number of cross-links per chain is kept 
fixed to distribute the cross-linkings over all chains. The interaction between the polymeric beads and the 
liquid particle follows the modified Lennard-Jones potential with the interaction parameter called a p . 

In order to make a two dimensional setup, we fix the size of the simulation box in z direction. The 
minimum depth of the simulation box is determined by two considerations. Lennard-Jones particles should 
not interact with themselves through the periodic boundaries. This condition is met by choosing the depth 
to be twice the cutoff radius. Since the cutoff radius is 5 d, 10 d is sufficient to remove the Lennard-Jones 
self-interaction of particles. Second, polymeric chains should not interact with themselves through periodic 
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boundaries. In our simulations, polymer chains consist of N = 32 monomers which makes their average 
radius of gyration equal to R g = yj^ 1 = 2.3 d. Therefore 2 (R g + r c ) = 14.6 d removes most of the self¬ 
interaction of the polymer chains. This condition is only important for polymeric liquids which have a small 
number of cross-links since the forces due to FENE cross-links are much stronger than the Lennard-Jones 
interactions. Thus the dynamics of the gel is mainly determined by the FENE-links. In this study, the 
depth of the setup is chosen to be 14.7 d. In addition, there is a maximum value for the depth to avoid the 
capillary instability. For a radius of R sa 75 d, capillary instability happens when depth of the simulation 
box is approximately 475d which is far above the considered range for the depth. 

In order to create the polymeric particle, first a big chunk of polymeric chains is made. Then those chains 
are positioned inside a bath of solvent, which naturally leads to the formation of a (cylindrical) polymeric 
droplet. Based on the desired cross-linking density, the polymeric droplet is cross-linked using (131 and then 
relaxed to reach its final equilibrium size. In our simulations, the resulting particle has 250,000 beads and 
its radius is 75.3 d ± 0.3d. 


C. Modelling of interface and particle 

After preparation of the equilibrated interface and cross-linked particle, the particle is inserted at the 
interface by carving out some volume of the liquid. It has been observed that a gap could be produced 
where the three phases are in contact which can be removed by increasing the pressure normal to the 
interfac^ 4 ^. After the insertion, the system is equilibrated under the NVT ensemble for 20, OOOr. The Nose- 
Hoover thermostat is used to keep the temperature constant at 300 K. We use the time step dt = O.Olr. After 
the equilibration, the shape of the particle is measured every 10, OOOr measurement steps. The criterion for 
convergence is that the shape of the particle in the three consecutive measurement steps is identical. 



FIG. 4 Density contours ip, as defined in (141. To obtain ip contours, the density is computed inside the square 
boxes of size 0.9d 2 in the x-y plane as depicted by the white grid. In addition, for the measurement of the 
macroscopic contact angle 9 S , a region of size 5 d is excluded from the contact point to exclude the effect of finite 
width of the fluid interface (shown with the line contours). 


The main goal of the study is to quantify the shape of the particles as a function of the various parameters. 
To obtain the shape from a molecular dynamics simulation, first the simulation data is extracted every 20r 
and the simulation box is shifted in a way that the center of mass of the particle remains fixed in all output 
frames. Then the simulation box is divided into smaller square boxes with an area of (0.9d) 2 in the x-y 
plane, as it is depicted by the white grid in Fig. [4] Inside each box the densities of polymeric chains and 
two liquid phases are calculated for each frame. Next, the densities are averaged over 500 frames. Finally, 
iso density contours for each measurement step are obtained usingpH 

, Ps - {pi + P 2 ) 

ip = -, 

Pi + Pa + P2 


( 14 ) 
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where p a and pb denote the density of two liquid phases a and b and p s denotes the density of the particle. 
A typical density contour is shown in Fig[4j Particles have four symmetric quadrants which can be used to 
increase the amount of data available for the averaging process by a factor of four. To obtain a macroscopic 
contact angle from the simulations, a small region on the order of the interface thickness, here chosen 5d, 
needs to be excluded. A typical contact angle measurement is shown in Fig. [4] 


D. Calibration 

In order to compare the results of the simulations with any macroscopic model, it is required to calculate 
the Young’s modulus E of the particle and the surface tensions q s and 75 . This calibration is outlined below. 


1. Elastic modulus measurement. 

To measure the Young’s modulus we prepare a cubic box filled with polymers with the same cross-linking 
density as the particle. Then, the pressure of the cross-linked polymer box is equilibrated independently 
in x : y and z directions in order to relax stresses to zero in all directions. Then the box is stretched in one 
direction by imposing a very small amount of strain, typically e s < 0.02. The length of the box is kept fixed 
in the stretched condition and the stress a s in that direction is measured. The lateral stresses are kept equal 
to zero during the stretching and relaxation stages, so that one directly probes the Young’s modulus as 


E m = —, (15) 

e s 

where At shows the time interval over which stress is averaged. Fig. [5] presents the normalized value 
of the Young’s modulus for different intervals of averaging. In molecular dynamics, stress exhibits much 
larger fluctuations as compared to velocity and density, since it is directly obtained from integrating rapidly 
changing interatomic forces (see e.g. Razavi et alpd). In our simulations, the final value of elasticity is 
obtained by averaging the values of elasticities for 50,000r to 100, OOOr measurement periods, where we 
note that the averaging needs to be longer as the stiffness of the gel decreases. 



At 


FIG. 5 Variation of the instantaneous Young’s modulus E&t (normalised by its time-averaged value E) with the 
interval of averaging, for two systems of different E. Error bars show the order of the fluctuations for each interval 
of the averaging, and give an indication of the convergence. For smaller average elasticity, a longer period of 
averaging is required to measure the Young’s modulus. 
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2. Surface tension measurement. 


The surface tension of the particle-liquid interface is measured using the Kirkwood-Buff formula for a 
planar interfaces. 

This is done by creating a seperate setup with a planar interface between the liquid phase and the cross- 
linked polymer chains. Fig. [ 6 ] shows how the resulting surface tensions vary with cross-linking density. At 
small cross-linking densities (p c < 1 ), the surface tension increases weakly with the cross-linking density and 
its value is very close to the liquid droplet case which has no cross-linkings. As p c increases, its effect on the 
surface tension becomes more significant. 



FIG. 6 Variation of the surface tension of the gel-liquid interface, 7 s , with the density of the cross-linking for two 
sets of interaction parameters: blue symbols are for a = 0.5, and a p = 0.4 and black symbols are for a = 0.8, and 
a p = 0.8. The planar interface is prepared inside a 120d x 120d x 120d box. After 10000 t of equilibration, the 
surface tension is measured during a 20000 r measurement period. 


IV. RESULTS 

A. Partially wetting particles 


We first consider simulations in the partially wetting regime, for which 7 > 1/2. In addition, we focus on 
situations where the contact angle is fairly large ( 9 S > 70°), so that deformations are small and fall within 
the linear regime. The key result is shown in Fig. [TJ where we characterise the particle shape by the radius 
r p as a function of the polar angle <f>, for a softness varying from S = 0.04 (fairly rigid, high cross-linking 
density) to S = 00 (liquid). Theoretical curves are shown as the solid lines and the results of the molecular 
dynamics simulations are shown using symbols. Except for the droplet (S = 00 ) which has a spherical cap 
shape, theoretical curves are computed directly from Eq. 0. The agreement between the macroscopic 
theory and simulation is excellent, and emphasises the strength of the continuum approach even for discrete 
and fluctuating systems at small scales. 

A few observations are of interest. For S = 0.04, the deformation with respect to an undeformed particle 
is restricted to a small region close to the contact line, up to about 4> ~ 5°. The range of the deformation 
increases with S, and extends to the entire droplet at large S. Also the magnitude of the deformation 
increases, and as a consequence of the large deformations one observes some small deviations from the linear 
theory (e.g. for the blue squares). In the limit of very large softness, however, the elasticity drops out of 
the problem altogether and one finds again a perfect agreement with the purely capillary shape of a liquid 
droplet. 
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FIG. 7 Shape of the particle in polar coordinates (r, <j >). Symbols show the results of the molecular dynamics 
simulations. Solid lines are the theoretical predictions obtained from Eq. |7|. For S = oo, a liquid droplet at a fluid 
interface, the line corresponds to a spherical cap shape. 


Next we investigate the maximum deflection A R in order to further quantify the particle deformation. An 
important prediction of the macroscopic theory is that the scaled deflection AR/ ARd rop is only a function of 


the softness S, and not of 7 . The prediction (111 is compared directly to the molecular dynamics simulation 
in Fig. [ 8 j as solid line and symbols respectively. Again, a perfect agreement is found between the continuum 
theory and molecular simulations. At low S the deformation radius increases ~ S'log (S'), as shown by (10), 
while at large S the deflection saturates to the value known for a purely liquid drop AR^ rop . As expected, 
the crossover between the regimes of small and large deformations arises when S ~ 1. 



10 2 10 3 


FIG. 8 Maximum radius of an elastic particle A R (normalised by that of a liquid drop A Rd rop ), as function of 
softness S. Symbols correspond to molecular dynamics simulations, the solid line corresponds to the prediction 
given by The far right data point corresponds to a liquid droplet (S = 00 ). The radial extension increases 
continuously with S, and saturates at the value of a liquid droplet. 


Finally, we can verify the prediction of the macroscopic analysis in section [TT] that the contact angle is 
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determined by the surface tensions according to the Neumann balance, regardless of the Young’s modulus 
of the particle. In Fig. [9j we plot the contact angle measured in our simulations 0 S rescaled by the Neumann 
angle 0<j rop . The ratio is indeed very close to unity, over the entire range of softness - even for the stiffest 
particles - meaning that the contact angle of the solid is governed by the surface tensions. Data in Fig. [9] 
covers contact angles from 40 to 82 degrees, meaning that the Neumann balance is valid even for the contact 
angles much smaller than 90°. 



FIG. 9 Variation of the normalized contact angle 9 a /8dro P with the softness. The solid line shows the theoretical 
value ( 6 S = 8drop) and symbols are the results from the molecular dynamics simulations for 8 S between 40° and 82°. 
The simulation results fall within the 4% of the theoretical values. The scatter in the data is because the contact 
angle does not correlate with the softness. 


B. Completely wetting particles 


For the partially wetting regime we have seen that the pulling force of the fluid interface is resisted by 
both the surface tension of the particle and its bulk elasticity: the capillary balance determines the geometry 
near the contact line, while at distances larger than 7 a /E the bulk elasticity is predominant. An interesting 
paradox arises when the surface tension of the solid is not sufficiently strong to balance that of the fluid- 
fluid interface. For the uncross-linked polymers this leads to a wetting transition, leaving a thin film at the 
interface, but this is clearly prohibited when the polymer chains are strongly cross-linked. As we will argue 
in section [Vj the situation of the complete wetting appears to be the generic case for microgel particles. In 
Fig. 10 we present two cases with equal Young’s modulus and ratio of surface tensions, i.e. identical S and 
7 , but the bath surface tension is so high that both particles are in the complete wetting regime (7 < 1 / 2 ). 
Both particles also have an identical homogeneous density of cross-links, but the difference between panel 
(a) and (b) is the topology of the cross-linking network: in the upper panel all polymeric chains are strongly 
entangled, while lower panel there are free chains in the particle that are not cross-linked. Note that even 
the spatial distribution of the links is similar between the two cases. The snapshot in Fig. [To| r shows that 
free chains get pulled out of the particle and produce a polymer film at the interface. The particle without 
the free chains, by contrast, is very strongly stretched (Fig. 10 1 ) - yet due to the finite extensibility of 
the FENE-polymers, the particle remains at a finite size. In this case, the excess force near the contact 
line 7 b — 27 s is balanced by the strong nonlinear elasticity of the cross-linked gel. Furthermore, there is an 


asymmetry between the left and right hand side of the particle in Fig. 10 a) which also could be observed for 
the soft particle in complete wetting regime in Fig. |T| This shows that for particles in the complete wetting 
regime, minor asymmetry in the cross-linking can affect the symmetry of the particle, highlighting the role 
of microscopic properties on the shape of the particle, while for particles in the partially wetting regime, 
surface tension undermines such effects. This result illustrates that (i) the completely wetting regime cannot 
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be captured by linear elasticity, and (ii) contrarily to the partially wetting case, molecular details become 
of key importance to describe the equilibrium shape. 



FIG. 10 Effect of the cross-linking topology on the deformation of a particle in the complete wetting regime, (a) 
Snapshot of a particle without any free chains. The particle stretches until it reaches a final size, (b) Snapshot of a 
particle with 30% free chains. Free chains are pulled out of the particle producing a thin liquid film. Both particles 
have identical thermodynamic properties, S = 0.45 and 7 = 0.37. 


Finally, let us investigate the influence of the softness on the particle shape in the case where no free 
chains are present. This is revealed in Fig. 11, where we show the relative extension of the particle A R/R 
as a function of S. The snapshots at small S reveal how a thin region near the contact line is pulled out of 
the particle: due to the complete wetting condition, the pulling force of the bath surface tension cannot be 
balanced by the surface tension of the gel, and hence no Neumann triangle can be formed. As a consequence, 
the particle elasticity must provide a highly localised stress at the contact line, which is achieved here be 
forming a cusp-like tip. Clearly, this feature is beyond linear elasticity and the linear theory does not suffice 
to capture the complete wetting regime. Our simulations show that the size of this tip-region increases as 
the particle gets softer, and in the range where S ~ 1 it invades the entire particle. In all cases, the strongly 
nonlinear deformation provides an elastic force that can compete with the surface tension force near the 
contact line. This is to be contrasted with the partially wetting case where elasticity does not play a role at 
the contact line. 


V. DISCUSSION 

We have investigated the shapes of deformable, elastic particles when adsorbed at a fluid-fluid interface. 
Apart from the elasticity of the particle, we have demonstrated that one crucially needs to consider the 
surface tension of the particle. In particular, one needs to distinguish cases of partial wetting and complete 
wetting. It was found that partially wetting particles simulated by Molecular Dynamics are accurately 
described by the continuum elasticity theory, while for the case of complete wetting, molecular details such 
the cross-linking topology turn out to be important. While the present paper has focused on symmetric 
particles, with equal affinity to the two liquid phases, it is interesting to consider more general values for 
the surface tensions. The case of partial wetting is reached when the sum of surface tensions of the solid 
7 s i and 7 S 2 , respectively with fluid “ 1 ” and “ 2 ”, are larger than that of the bath 7 ^. At the same time, for 
adsorption to be energetically favourable with respect to the particle being in phase “ 2 ”, one requires 7^1 
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FIG. 11 Maximum radius AR (normalised by the initial drop radius R), as a function of S in the case of complete 
wetting. All polymers are cross-linked so that the particle reaches an equilibrium with finite deformation. All 
particles are subjected to the same pulling force and are in the complete wetting regime. For numerical reasons, the 
largest value of S was achieved by a slightly smaller particle. 


to be lower than 75 + 7 s2 . Partially wetting and adsorption is thus reached when the following inequality is 
satisfied: 


7b - 7s2 < 7s 1 < 7b + 7 s2- 


(16) 


It is interesting to consider this condition in the context of microgel particles. These are hydrogels that 
have a nonzero but small surface tension with respect to the water phase in which they are immersed. 
Denoting the water as phase “2”, we typically find that the gel-water surface tension 7 s2 is much smaller 


than the other surface tensions. In terms of the inequality (16) this implies that the window for the partial 
wetting regime of adsorbed particles is extremely small, and that one generically expects the particle to be 
in completely wetting situations similar to what can be seen in our Fig. |10| Exper imentally, it has remai ned 
challenging to explain the core-corona like shape found for the microgels 18 23 31 . As is shown in Fig. 10 


the core-corona shape could happen either due to the deformation of the entangled polymeric network of 
the microgel particles or due to the adsorption of the individual polymeric chains. For a typical microgel 
particle, with a Young’s modulus of 100 KPa' 1 , a particle-oil surface tension of 50 mN/m, and a typical 
radius of a 0.5 micron, the softness parameter will be of order unity. Given that experiments are likely to 
be in the complete wetting regime, according to Fig| ll[ th e expected deformation is therefore of the order 
of the particle radius consistent with the experiments^® 31 . The important observation from comparing the 
snapshots of FigJTT] is that for a particle without free chains, the size of the corona (cuspy shape close to the 
contact line) which is clearly present for small S (snapshots a to c) vanishes for the particle whose softness is 
close to one (snapshot d ) which is the typical experimental value. This shows that for S ss 1, the deformation 
spans over the entire particle shape instead of producing an elongated cusp close to the contact line. If we 
compare this observation to the experimental observations that the corona size is comparable to the particle 
core siz d 26 * 31 !, we can conclude that the details of adsorption of individual polymer chains play an important 
role in producing a core-corona shape for microgels at a fluid interface. It should be noted that in this study 
as the first step to model the behavior of microgel particles at a fluid interface, particles are considered to 
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be in the dry state which means that the solvent particles do not penetrate into the microgels and the effect 
of particle swelling is ignored. It would be interesting to extend the present work to incorporate the effect 
of particle swelling which is usually present in experiments. 


VI. CONCLUSIONS 

The shape of a soft particle at a fluid interface is studied both numerically and analytically. It is shown 
that the surface tension of the particle, expressed in terms of the dimensionless parameter 7 plays an 
important role in characterizing the particle shape. For 7 > 1.4, i.e. particles whose surface tension is 
large compared to the fluid surface tension, it is analytically shown that the contact angle of the particle is 
equal to the Neumann contact angle of a liquid droplet with similar 7 . In addition, an analytical formula is 
derived for the particle shape and it is shown that it scales with the shape of the liquid droplet with similar 
7 where the scaling coefficient is only a function of the particle’s softness S = The derived analytical 
equation for the particle shape and contact angle is examined using the coarse grained molecular dynamics 
simulations of particles made out of entangled polymeric chains and it is shown that the analytical and 
numerical results are in a perfect agreement. For 0.5 < 7 < 1.4, molecular dynamics simulations show that 
the analytical formula for the particle shape loses its accuracy as it is expected but the contact angle of 
the particle is still equal to the Neumann angle. For 7 < 0.5, which considers the particles in the complete 
wetting regime, it is also shown that the shape of the particle depends on the microscopic details of the 
cross-linking. For example, two particles with similar macroscopic properties but different number of free 
chains produce different final shapes. Using the results of the particle in the complete wetting regime, it is 
argued that the core-corona structure which is observed in the majority of the microgel experiments appears 
due to the adsorption of free chains inside the particle to the interface. 
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Appendix A: Appendix 

As we are dealing with a two-dimensional problem, we can exploit the Airy stress function <I> and solve 
the biharmonic equation V 4 4> = (P. The corresponding stress and displacement fields can be determined 
as partial derivatives of $. In polar coordinates, we can proceed by a Fourier expansion using the classical 
solution by Michell. The up-down symmetry is such that we only need to consider Fourier modes of the 
type cos ruj), where n = 0, 2,4, •••. The isotropic solution n = 0 comes with powers r 2 , r 2 lnr, lnr. The 
n = 2 ,4, 6 , • • • solutions to the biharmonic equation come with radial powers of the type r n , r n+2 , r ~ n , r _n+2 . 
Demanding regularity of stress and strain at r = 0, we cannot allow for r~ n , r _n+2 or the logarithmic terms, 
so that: 


4>(r, (j>) = Ar 2 + ^ cos n</> ( a n r n + b n r n+2 ) , 

n =2 

The corresponding displacements (taken at the free surface r = 1) become 57 : 


u r {l,4>) 

^( 1 , 0 ) 


3 S 

27 

3 S 
27 


OO 

cos rup (—na n — nb n ), 

n =2 
00 

y sin ruj) ( na n + (n + 2 )b n ), 

n—2 


(Al) 


(A2) 

(A3) 
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where S is the anticipated dimensionless parameter to characterise the softness and 7 is the ratio of the sur¬ 
face tensions. Due to the incompressibility (v = 1/2) there is no isotropic contribution to the displacement. 
Likewise, one expresses the stress at r = 1 from the expansion (Al I as 


OO 


<r rr (l, (f>) = 2 A + cos n4> [(n — n 2 )a n 

n=2 


+(n + 2 — n 2 )b n ] , 

(A4) 

cr r 0 (l, 0 ) = Ynsmmf, [{n - l)a n + (n + 1 )b n ]. 

(A5) 


n—2 


The final step is to expand the boundary conditions ( 5j6 ) in a Fourier series and compare this to the 
representation ( A4|A5 ). This will give two equations for the coefficients a n ,b n . In particular, the normal 
stress boundary condition becomes, using also (A2| 


OO 

0 >( 1 ,4>) = — + y; c n cos w</>, 

n—1 


(A 6 ) 


with coefficients 


2 2 3 

Co = — , Cn = - + -S [n(n 2 - 1 )(a„ + b n ) 1 n = 2,4,6, •• • , 

7T 7T 2 

while c n = 0 for odd n. Solving the coefficients a n , b n finally gives 


" tt (1 - n) [1 + | Sn] 

The displacement field can be summarized in explicit form: 


’ " 7r(l + n) [l + ^Sn] 


. ,. 3S v/k oI , 

AA <P) = 2^cos2fc</ 


r 2fc+i + 


„2fc-l 


2 k 


2 k 


(1 - 2k) (1 + 3Sk) 


3 S 


H(r, <t>) = g— 51 sin 

9 ,?r fc=i 


(l + 2fc) (1 + 3S7c) y ' 
2 k 


Jlk— 1 


(1-2*0(1 + 35*) 


2fc+l 


2k+ 2 


(1 + 2fc) (1 + 3Sk) 
The shape of the deformed particle is thus given by 


r P (</>) = l + u r (l,<f>) 

3 S y, 4kcos2k(j) 

Or-j'Tr ‘ f 1 — . 


2y7r (1 — 4fc 2 )(l + SS'fc) 


(A7) 


(AS) 


(A9) 


(A10) 


(All) 
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